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ABSTRACT 

Together with the variational indicators of chaos, the spectral analysis methods have 
also achieved great popularity in the field of chaos detection. The former are based 
on the concept of local exponential divergence. The latter are based on the numerical 
analysis of some particular quantities of a single orbit, e.g. its frequency. In spite of 
having totally different conceptual bases, they are used for the very same goals such as, 
for instance, separating the chaotic and the regular component. In fact, we show herein 
that the variational indicators serve to distinguish both components of a Hamiltonian 
system in a more reliable fashion than a spectral analysis method does. 

We study two start spaces for different energy levels of a self-consistent triaxial 
stellar dynamical model by means of some selected variational indicators and a spectral 
analysis method. In order to select the appropriate tools for this paper, we extend 
previous studies where we make a comparison of several variational indicators on 
different scenarios. Herein, we compare the Average Power Law Exponent (APLE) and 
an alternative quantity given by the Mean Exponential Growth factor of Neary Orbits 
(MEGNO): the MEGNO's Slope Estimation of the largest Lyapunov Characteristic 
Exponent (SEZLCE). The spectral analysis method selected for the investigation is the 
Frequency Modified Fourier Transform (FMFT). 

Besides a comparative study of the APLE, the Fast Lyapunov Indicator (FLI), 
the Orthogonal Fast Lyapunov Indicator (OFLI) and the MEGNO/SEZLCE, we show 
that the SEILCE could be an appropriate alternative to the MEGNO when studying 
large samples of initial conditions. The SEfl^CE separates the chaotic and the regular 
components reliably and identifies the different levels of chaoticity. We show that the 
FMFT is not as reliable as the SEfl^CE to describe clearly the chaotic domains in the 
experiments. We use the latter indicator as the main variational indicator to analyse 
the phase space portraits of the model under study. 

Key words: methods: numerical - galaxies: kinematics and dynamics. 



1 INTRODUCTION 

The understanding of a realistic dynamical model is strongly 
related to the capability of identifying the chaotic and the 
regular nature of its orbits. 

One the one hand, one of the main features of chaotic 
orbits is their high sensitivity to the initial conditions (here- 
inafter i.e.), and the concept of local exponential divergence 
applies perfectly to identify such correspondence. On the 
other hand, the regular orbits do not show this behaviour. 
A seminal contribution to the field of chaos detection has 
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been made by Lyapunov when he introduce d the Lyapunov 
Characteristic E xponents, hereinafter LCEs (|Oselededl 19681 ; 
lLvapunov|[l99^ . The LCEs are theoretical quantities that 
measure the local rate of exponential divergence. Therefore, 
the LCEs are appropriate to distinguish between chaotic and 
regular motion. 

The LCEs constitute the backbone of variational 
indicators of chaos (VICs). Moreover , their numeri- 
cal implementation dBenettin et al.l Il980l ; iFroeschiej 1 19841 : 
iTancredi et al. | l200ll ; ISkokosI I2010I L i.e. the Lyapunov In- 
dicators (Lis), was a major contribution to the devel- 
opment of both fast and easy-to-compute new VICs. 
Nowadays, there is a plethora of VICs in the litera- 
ture, such as the Smaller Alignment Index, SALI jS kokos 
l200ll ; ISkokos et al]|2004 ISzell et alfiool iBountis fc Skokosl 
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l200d : ICarpinterd |2008l : lAntonopoulos et all l201Ch . and 
its gene ralized version, the Generalized Alignment Index, 
GAL I |Skokos et all 120071 . |200S| ; |m anos fc Athanassoulal 
l201l[ ); the Mean Exponential Growth fa c tor of Nearby 
Orbits, MEGNO (ICincotta fc Siinol l200d; ICincotta et "all 



20031; iGiordano fc Cincottal 12004 iGozdziewski et al.ll2005l; 



Gavon fc Boisl 20081; iLemaitre et alj 



201(1 ICompere et alj I2OIII ; IMaffione et al 



2009 



Fast Lyapunov I n dicat or , FLI 
Froeschle fc Legal 1 19981, 1200(1: 



Guzzo et al 



.Hince et al.l 
201 lal); 



the 



Froeschle ct al 
Le ga fc Froeschle 200.1 



1997allb 



20021: IFroeschle fc Legal l200d; iPaleari et all 



20081 : iTodorovic et al l I2OO8I ; iLega et all |2010| ); its 



ant, the Orthogonal Fast Lyapunov Indicator, OFLI 
|Fouchard et all |2002| ) : and a second order variant of t he 
OFL I, the OFLl| T (jBarrid 120051 : iBarrio et all l2009al lbl. 
l20ld ). Furthermore, we can include the VICs based on the 
properties of dynamical spectra like the invariant spectra of 
Stretching Nu mbers or Local Lyapunov Characteristic Num- 
bers, LLCNs (|Contopoulos et all 1 19971 ; ILega fc Froeschle! 
1998), the invariant spectra of H elicity or Twist An- 
gles dVoglis fc Contopoulosl 1994 ; IContopoulos fc Voglisl 

1997 ; IFroeschle fc Legall 19981; 



199d.ll997l;IContopoulos et al 



Lega fc Froeschlelil99Sl; IVoglis et alTll998h and the Spectral 
Distance ( Voglis et al.l 1 19991 ). Fi nally, we include the Rel- 
ative Lyapunov Indica tor, RLI (|Sandor et all l200d , |2004| . 
120071 : ISzell et alj|2004h . which is not based on the evolu- 
tion of the solution of the first variational equations as the 
others, but on the evolution of two different but very close 
orbits. 

Notwithstanding the large number of VICs cited, 
we could also include alternative methods, for in- 
stance, the Average Power Law Exponent, APLE 
jLukes-Gerakopoulos et all [2008), which is based on 
the concept of Ts allis entropy; evolutionary algorithms 
l|Petalas et alJ l2009): and the record could indeed continue. 

Together with the VICs, the spectral analysis meth- 
ods (SAMs), based on the numerical analysis of some 
particular quantities of a single orbit, e.g. its fre- 
quency, have also achieved great popularity in the 
field of chaos detection. Amon g the SAMs we find 
the method outlined by Laskar] dl990h: t he Frequenc; 
Map Analysis, FMA (jLaskar et al 



Map 
IPapai 



^tUI): t ne .frequency 
19921 ; lLaskarl 1 1993 ; 



iPapaphilippou fc La skar 1996, 1 19981 1 and a modification of 
the latter in order to improve the precision in the com- 
puted frequencies and amplitud es, the Frequency Modi- 
fied Fourier Transform, FMFT (jSidlichovskv fc Nesvornvl 
1997) or other alternatives using Fast Fourier Trans- 



form (FFT) technic 


mes ( 


Michtchenko & Ferraz-Melldll99d; 


Ferraz Mcllo ct al. 


2005 


: iMichtchenko et al.ll20ld). among 



others. As an easy identification of the important resonances 
(regions of phase space occupied by orbits with commensu- 
rable orbital frequencies) is possible with the SAMs, they 
are able to identify the major resonant orbit families. The 
relative importance of each family regarding the phase space 
can be assessed from the number of orbits associated with 
the particular resonance the family belongs to. That is, the 
SAMs are not only suitable to determine the resonance web 
like the VICs llKaneko fc Konishill 19941; ICincotta et al.l|2003l ; 



IFroeschle et all 200d ; Lukes-Gerakopoulos et al. 12008 ) but 
also to identify the individual resonances by their frequency 
vectors. 

Though the conceptual bases of the VICs and the 



SAMs are totally different, many researchers use both 
types of techniques for similar goals. Nevertheless, we 
show herein that this seems to be a no reliable approach 
to study some important aspects of the dynamics of a 
Hamiltonian system. Moreover, in the literature there are 
plenty of evidence that the complementary use of both 
types of techniques results an efficient and solid way to 
;ather dynamical information from a Hamiltonian system 
Froeschle et al.lll997bl: IFroeschle fc Legalll998l:lGuzzo et al.l 
20021 ; IKalapotharakos fc Voglisl |2005l; IBarrio et all l2009bK 
For instance, in IFroeschle et al.1 |l997b) the authors com- 
pare the sensitivity of the FLI and th e FMA on the stan- 
dard map and on the iHenon fc Heilesl (|l964l ) potential. In 
iGuzzo et al.l (|2002l ) the authors use the F LI and the Analyt- 
ically Filtered Fourier Analysis, AFFA (|Guzzo fc Benettir] 
120011) . another example of SAM, to determine the value of 
the critical parameter at which the transition from the unfor- 
tunately called Nekhoroshev to the Chirikov regime occurs 
in a quasi-integrable Hamiltonian model and a s tanda rd 4- 
dimensional map. In IKalapotharakos fc Voglisl (|2005l ). the 
authors use a measure related w i th the LCEs, the alignment 
indices l| Voglis et al.ll 19981 . Il999l ; ISkokosll200ll ) and an index 
computed with the FMFT to study two different N-body 
models simulating elliptical galaxies. 

One of the aspects that makes a research accurate is 
the reliability of the tools available for the analysis. This 
reliability is a measure of, for instance, the confidence and 
efficiency of such tools. Regarding chaos detection, we be- 
lieve that it is on the diversity of the employed techniques 
that the fruitfulness of the research relies. Therefore, we 
started a series composed of four papers in order to achieve 
the best combination of t ools to study a given d ynamical 
problem. In a first paper, IMaffione et al.1 (|2011al ). we test 
the MEGNO against the most widely used chaos detection 
technique in the literature, which is the LI, on a somewhat 
realistic proble m: a self-consistent triaxial stellar (ScTS ) dy- 
namical model (jMuzzio et al.l|2005l: ICincotta et a"il200Sl ). In 
the second and t hird papers. | Maffione et al.l (|201 lbT ) (here- 
inafter Mil) and lDarriba et al] (|2012l ) (hereinafter D12), we 
make a comparison of most of the VICs mentioned earlier on 
ra ther simple prob l ems, i.e. on mappings and on the model 
of|H enon fc Heilesl (|l964l ). respectively. Our purpose here is 
to extend the previous comparisons of VICs to alternative 
techniques, including a spectral analysis method. Further, 
this last paper of the series is regarded as a closing report 
on the subject and the conclusions comprise the whole in- 
vestigation. We continue with the ScTS model and study 
two start spaces for four different energy levels by means of 
several chaos detection tools. In order to select appropriate 
VICs for the study, we consider the results of Mil and D12 
and add a rather new technique, the APLE, and an alter- 
native quantity given by the MEGNO: the MEGNO's Slope 
Estimation of the largest LCE (SEflLCE) to their compari- 
son. The MEGNO is part of the suggested VICs to compose 
the CIsEQ in D12. In fact, we are interested in studying if its 
variant, i.e. the SEflLCE, is more reliable than the MEGNO 
on its own. 

The study is organized as follows: in Section [2] we de- 



1 CIsF: "Chaos Indicators Function" . It is a minimal and efficient 
package of variational indicators to study a general Hamiltonian. 
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Table 1. List of the most important acronyms in alphabetical 
order. 



Acronym Definition 



APLE Average Power Law Exponent 

CIsF Chaos Indicators Function 

FLI Fast Lyapunov Indicator 

FMA Frequency Map Analysis 

FMFT Frequency Modified Fourier Transform 

LCEs Lyapunov Characteristic Exponents 

Lis Lyapunov Indicators 
MEGNO Mean Exponential Growth factor of Nearby Orbits 

OFLI Orthogonal Fast Lyapunov Indicator 

SEiLCE Slope Estimation of the largest LCE 

ScTS model Self-consistent Triaxial Stellar model 

VICs Variational Indicators of Chaos 



Table 2. Adopted values for the coefficients of the functions /, 
given by Eq. Q. 







a 


a 


<5 


ac 


n 


= 


0.92012657 


1.15 


0.1340 


1.03766579 


n 


= X 


0.08526504 


0.97 


0.1283 


4.61571581 


n 


= z 


-0.05871011 


1.05 


0.1239 


4.42030943 



The ScTS model is an A^-body potential which is nei- 
ther smooth nor stationary. As the authors needed to com- 
pute the orbits and their corresponding Lis, they froze the 
potential and represented it with a quadrupolar approxima- 
tion. The equation that reproduces the potential is: 

V(x) = -/o(x) - /„(x) ■ (x 2 - y 2 ) - /,(x) • (z 2 - y 2 ) 

with x = (x, y, z) and where 

/»(*) = ^-BT . (!) 

[Pn n + Sn* 1 ] a " 

with p\ the square of the softened radius given by pi — 
x 2 + y 2 + z 2 + e 2 when n = 0, or p 2 n = x 2 + y 2 + z 2 + 2 ■ e 2 
for n — x, z, and a n , 8 n , a n , ac n are constants. 

The selected value for the softening parameter is e ~ 
0.01 for any n. The adopted values for the constants a n , S n , 
a n and ac n are given in Tabled 

The stationary character of the parameters given in Ta- 
ble[2]were tested by performing several fits at different times 
after virialization, resulting in a preci sion of 0.1%. For fur - 
ther details on t he Sc TS model refer to lMuzzio et~aLl (|2005h : 
ICincotta etafl (|2008i ). 

We have already used t he same model to t est the 
MEGNO in a previous study (jMaffione et alj|2011ef ). 

The spaces of i.e. used in this investigation are the 
Px a — Pz and the xo — zo start spaces for four energy sur- 
faces: —0.1, —0.3, —0.5 and —0.7. The p Xo — p ZQ start space 
is a set of i.e. of the form (p X0 ,p ZQ ) with xo = yo = 
zp = while p vn is restraine d by the energy condition 
|Papaphilippou fc La skar 199c]). The xo — zq start space is 
a set of i.e. of the form (xo, zo) with yo = 0, p xn = p zn = 
and p vn is restrained by the energy condition (|Schwarzschildl 
0993 ). On the one hand, the combination of both start spaces 
is an adequate choice to sample many of the different or- 
bital families within the triaxial model l|Schwarzschildlll993l : 
IPapaphilippou fc La skar 1998). On the other hand, the sta- 
ble boxlets (i.e. resonant box orbits) which are centrophobic 
are not well represented in our p XQ — p ZQ start space. 

Therefore, the p XQ — p ZQ and the xo — zo start spaces of 
the ScTS model seem to provide fairly realistic scenarios for 
testing many characteristics of both the VICs of the package 
and the FMFT. 

3 THE TECHNIQUES 

The aim of this section is to briefly introduce the techniques 
we will use in the forthcoming study. 

3.1 The LI, the MEGNO and the SEZLCE 

Considering a continuous dynamical system defined on a 
differentiable manifold S, where <& (x) = x(f) characterises 



scribe the ScTS model and introduce the p X(l — p ZQ and the 
Xo — zo start spaces. In Section we briefly introduce the 
chaos detection techniques to be used within the reported re- 
search. We start with the VICs: the LI, the MEGNO and the 
SEZLCE (Section EU); the FLI and the OFLI (Section [oT^ 
and the APLE (Section O. We finish with a short FMFT 
introduction (Section 13. 4[) . In Section [H we deal with the 
capability of the FMFT to resolve the different structures 
in a divided phase space (Section 14. 1[) and compare its per- 
formance with the SEflLCE in both start spaces of the ScTS 
model (Section 14. 2p . In Section [S] we study the efficiency of 
several VICs on the p X(l — p ZQ start space to decide which 
ones are appropriate to investigate the dynamics of the ScTS 
model (Sections l5.1l and l5.2[) . In Section[S]we analyse for four 
different energy surfaces, the p xg —p ZQ (Section l6.1|) and the 
xo — zo (Section 16. 2p start spaces by means of a comple- 
mentary use of the VICs selected in the previous section. 
Finally, Section summarizes the main conclusions of the 
whole investigation. 

As the extensive used of acronyms in the text might be 
overwhelming for the reader, in Table [1] we list in alphabet- 
ical order the most important acronyms used in this paper. 



2 THE MODEL AND THE START SPACES 

In order to obtai n a fairly realistic dy namical model of an el- 
liptical galaxy, in lMuzzio et al.l (|2005t ) the authors follow the 
cold disipationaless collapse of 100000 particles randomly 
distributed following a 1/r law within a sphere. The veloc- 
ities were randomly chosen from a spherical Gaussian dis- 
tribution but only their tangential component was retained. 
After the system had relaxed, there remained 86818 particles 
resembling an elliptical galaxy: the ScTS mo del. The sys- 
tem o bserves the de Vaucouleurs law shown bv lMuzzio et al.l 
(2005) in their Fig. 2. The model reproduces many dynami- 
cal characteristics of real elliptical galaxies, su ch as mass dis - 
tribution, flattening, triaxiality and rotation (|Muzzidl2006| y 
The ScTS model has a strong triaxiality and a flattening 
that increases from the border of the system to its center (see 
Table I in the same paper). The resulting triaxial potential 
has semi-axis X, Y, Z satisfying the condition X > Y > Z, 
and its minimum, which is close to —7, matches the origin. 
The potential is less flat tened than the mass distribution, as 
expected. See Table I in lMuzzio et alj (|2005| k 
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the state of the system at time t, being the state of the 
system at time t = 0, x(0) = xo. The state of the system 
after two consecutive time steps t and t' will be given by the 
composition law: <I> t+t =$'o$' . 

Moreover, the tangent space of x maps onto the tangent 
space of 3?'(x) according to the operator d x & 1 and following 
the rule w(f) = d x $ t (w(0)), where w(0) is an initial devia- 
tion vector (hereinafter i.d.v.). The action of such operator 
at consecutive time intervals satisfies the equation: 



If we suppose that our manifold S has some norm de- 
noted by || • ||, we can define the useful quantity: 



K**(w)|| 



called "growth factor" in the direction of w. 

Let H(p, q) with p, q £ be an iV-dimensional 
Hamiltonian, that we suppose autonomous just for the sake 
of simplicity. Introducing the following notation: 

x = (p,q) G K 2iv , f(x) = {-dn/dci, dH/dp) £ R 2N , 

the equations of motion can be written in a simple way as 

x = f(x). (2) 

Let 7(xo; i) be an arc of an orbit of the flow Q over a 
compact energy surface: Mh C R 2JV , Mh = {x : W(p, q) = 
h}, then 

7(xo;t) = {x(x ;t') : x £ Mh, 0^i'< t}. 

We can gain fundamental information about the Hamil- 
tonian flow in the neighborhood of any orbit 7 through the 
largest LCE: 1LCE, defined as: 



X~t = x[l(xo;t)] = lim -lnAt[7(x ;t)] 

t— >oo I 



with 

A t [7(x ;t)] 



(3) 



||d 7 *'(w)|| 



where ^^^'(w) || is an "infinitesimal displacement" from 7 
at time t. The fact that the 1LCE (and its truncated value, 
the so-called LI= linu_>T f In A t [7(xo; t)] for T finite) mea- 
sures the local mean exponential rate of divergence of nearby 
orbits is clearly understood when Eq. ([3| is written in an in- 
tegral fashion: 



X~t = lim 7 



IK*' (w) 
lld»#*'fw) 



dt' = 



|d 7 **'(w) 



(4) 



||d,**(w) 

where the bar denotes time average. 

Finally, the orbit is classified as regular or chaotic if the 
LI tends to zero (as ln(t)/t) or to a positive value, respec- 
tively (a detailed discussion on the theory and computat ion 
of the LCEs can be found in the review of lSkokod j201Ch ). 

Now, we can introduce the quantity Y [7(xo; t)] through 
the expression: 

which is related to the integral in Eq. @; i.e., in case of 
an exponential increase of ||ci 7 $ t (w)||, ||d 7 «I>'(w)|| = ||w|| • 
exp(x 7 i), the quantity Y 1 can be considered as a weighted 



variant of the integral in Eq. (Q. Notice that the quantity 
X 7 = Y-y/t verifies x-i ~ 2/t f° r regular motion and \i ~ X~i 
for chaotic motion, with t — > 00, which show that, in case 
of regular motion \i converges to faster than \i does 
(which goes to zero as lnt/t), while for chaotic motion both 
magnitudes approach the positive 1LCE at a rather similar 
rate. 

Instead of using the instantaneous rate of increase, Xi, 
we average the logarithm of the growth factor, ln(At) = \ 1 t. 
Introducing the time average (MEGNO): 



Y 7 = Y[ 7 (x ;*)] = i 



Kydt', 



we notice that the time evolution of the MEGNO can be 
briefly described in a suitable and unique expression for all 
kinds of motion. Indeed, its asymptotic behaviour can be 
summarized as follows: 



Y, 



(5) 



where a 7 = and b 1 w 2 for quasi-periodic motion, while 
a i = X7/2 and 6 7 m for irregular, stochastic motion. Devi- 
ations from the value 6 7 w 2 (and a 1 — 0) indicate that 7 is 
close to some particular objects in phase space, being 6 7 < 2 
or b 1 > 2 for stable periodic orbits (resonant elliptic tori), 
or unstable periodic orbits (hyperbol ic tori), respect i vely. 

Refer to Fig. 1 (c) and (d) from ICincotta etall (I2003T ) 
to see the general behaviour of the MEGNO for regular and 
chaotic orbits. 

It is possible to estimate the LI of the orbit from 
Eq. (0 if we apply a simple linear least-squares fitting o n 
the MEGNO l|Cincotta fc Simdl200d : fCincotta et al.ll2003l 1. 
This estimation is the earlier mentioned SEflLCE indicator. 
The least-squares fitting at the end of the integration pro- 
cess uses the last 80% of the points in order to avoid the 
initial transient. Thus, the SEflLCE almost considers the full 
history of the orbit which makes it a very sensitive VIC. 

The performances of the MEGNO's time evolution 
curves to characterise the different levels of stability and 
chaoticity of the orbits shown in Mil and D12 was excel- 
lent. On the other hand, the description of the regular com- 
ponents in a divided phase space by means of the final val- 
ues of the MEGNO (i.e. the values of the indicator at the 
end of the total integration time), was not satisfactory. The 
SE/LCE seems to be an appropriate alternative to solve the 
problem. 

Fig. [T] shows the performances of the LI, the MEGNO 
and the SEflLCE for two sets of i.e., one of them located 
in the p XQ — p ZQ start space (top panel) and the other set 
located in the xo — zo start space (bottom panel). In order 
to ensure that the VICs are well computed for a given or- 
bit, the total integration time should verify T >> T C (E), 
where the function T C (E) is some characteristic time-scale 
which depends on the energy surface E. The function T C (E) 
is taken as the period of the stable X-axis periodic orbit 
and we consider the energy surface E = —0.7 for both sets 
of i.e. (T c (— 0.7) ~ 7 units of time -hereinafter u.t.- see 
iMaffione etall |2011al )). Therefrom, on fixing the condition 
T > 10 3 T C (E) to obtain reliable values for the VICs, we 
conclude that a total integration ti me of 7 x 10 3 u.t. woul d 
be appropriate for the experiment (Maffi one et al.l [2011aT l. 
In the p XQ — p ZQ start space, the set is composed of 790 i.e. 
of the form p xo £ [0 : 2.5] and p ZQ = 0. In the xq — zq start 
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Figure 1. Final values of the LI, the MEGNO and the SEZLCE. 
The i.e. were taken in the p XQ — p ZQ start space (top panel) and 
in the xq — zq start space (bottom panel). Both sets of i.e. were 
taken on the energy surface —0.7 and the VICs were integrated 
for 7 X 10 3 u.t. The values of the indicators are in logarithmic 
scale. 



space, the set is composed of 646 i.e. of the form 20 £ [0 : 0.7] 
and xo = 0. 

In Fig.[T] the SEZLCE and the LI show very similar final 
values for chaotic orbits. The MEGNO grows exponentially 
fast and reaches higher final values (> 2) than the SEZLCE 
or the LI. The MEGNO clearly distinguishes between dif- 
ferent levels of chaoticity. Consequently, the SEZLCE seems 
to be an unnecessary alternative to the MEGNO for de- 
scribing the chaotic component. In case of a divided phase 
space, where regular and chaotic components are mixed, the 
asymptotically behaviour of MEGNO's final values for reg- 
ular orbits may hide their different stability levels (see Mil 
and D12). The SEZLCE tends to zero for regular orbits. 
Moreover, the SEZL CE tends to decrease fa s ter than the LI 
for r egular motion (|Cincotta fc Simicl |2000| ; ICincotta et al.l 
2003). So, no matter if the phase space is strongly divided, 
this behaviour allows the SEZLCE's final values to show cer- 
tainly (and faster than the LI) different levels of stability for 
the regular orbits. Thus, it seems to be appropriate to use 
the time evolution curves of the MEGNO for the individual 
analysis of the orbits (Mil and D12) and the final values 
of the SEZLCE as a global VIC to obtain the portraits of 
divided phase spaces (see Sections 14.21 16.11 and I6.2|l . 

The MEGNO has a low computational cost and the lin- 
ear least-squares fitting in the computation of the SEZLCE 
does not require a significant amount of extra time. Hence, 



the SEZLCE is an alternative quantity that has also a low 
computational cost. 

3.2 The FLI and the OFLI 

The F LI is a qua n tity i nti mately related t o the 
1LCE (IFroeschle et all Il997bllal IFroeschle fc Legal l200fj|; 
Lega fc FroeschlefeOoJlGuzzo et al.ll2002|; IFroeschle fc Lega ' 



20061 ) and the MEGNO jMestre et al.1 l201ll) . The FLI is 
able to distingu ish between regular a nd chaotic (weakly 
chaotic) motion (|Froeschle et al.lll997bl a[) and also between 
resonant and non-resonant motion (IFroeschle fc Legall2000l ; 
iLega fc Froeschlej I200J ; iGuzzo et all |2002| ) using only the 
first part of the computation of the 1LCE. 

Here, we use the ini tial definition of the FLI given 
in IFroeschle et al.l |l997bl lal) like in D12, i.e. we use a ba- 
sis of i.d.v. for i t s com putation. The definition given in 
IFroeschle fc Legal {2000) uses only one i.d.v., the results do 
not chan ge. Nevertheless, the co mputing time is obviously 
reduced ^Froeschle fc Legall2000l ). 

On an iV-dimensional Hamiltonian H, we follow the 
time evolution of 2N unit i.d.v. (e.g. the ScTS model is 3- 
dimensional and the basis consists of 6 i.d.v.). The FLI at 
time t is defined by the highest norm among the unit i.d.v. 
of the basis, as follows: 



FL/(t) = Bup[||w(t)i||,||w(t) a ||, 



||w(t) a 



For further details, refer to IFroeschle et~all (|l997bl lah. 

For both regular and chaotic motion, the FLI tends 
to infinity as time increases, but with completely different 
rates. The FLI grows linearly with time for regular motion 
and it grows exponentially fast for chaotic motion. Further- 
more, the FLI grows linearly with time, but with different 
rates in the case of resonant regular orbits and non-resonant 
ones. Although the curves corresponding to resonant and 
non resonant motions are separated, the oscillations due to 
the distortion of the orbits themselves may hide the distinc- 
tion between the l i bratio n islands and the tori. Therefore, in 
IFroeschle fc Legal l|2000l ) the authors not only introduce the 
FLI using only one i.d.v. but also replace the definition of 
the FLI by its running average. We also applied an average 
to smooth the curves of the FLI. 

In the case oftheOFLfl, we take the orthogonal compo- 
nent for the flow for each unit i.d.v. (w(t)f~ , i = 1, . . . , 2N) 
of the basis at every time step and retain the largest com- 
ponent among them: 

OFLI(t) = sup [w(t)i,w(t)a,---Mt)w] . 

This modification makes the OFLI a VIC which can 
easily distinguish periodicity among the orbits of the regular 
component. The OFLI for periodic orbits oscillates around 
a constant value, while for quasiperiodic motion and chaotic 
motion it has the same behavi our as the FLI. 

Refer to Figs. 5 and 7 from lFouchard et al.l (|2002l ) to see 
the general behaviour of the FLI and the OFLI for regular 
and chaotic orbits. 



2 Let us remark that the following definition of the OFLI 
is slightly different f rom the definition of the OFLI given in 



iFouchard et alj d2002t) where the authors use only one i.d.v. 
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3.3 The APLE 



In the early study of iTsallis et al.l l|l997h the authors show 
that when a nonlinear dynamical system is in the regime of 
the so-called edge of chaos, the q-entropy rate of increase 
remains constant for a long time (metastable states). More- 
over, they discuss that this behaviour may be associated 
with a power-law rather than with exponential sensitivity 
of the orbits to the i.e. The interested reader may refer to 
ITsallis et al. I (|2002h for a "pedagogical guided tour" on the 
subject . 

In iLukes-Gerakopoulos et all (|2008l ) the authors find 
that the time evolution of the i.d.v. for a weakly chaotic 
orbit 7(xo;t) (i.e. a combination of a linear and an expo- 
nential law: w(t) ~ c-yt + e yt , with c 7 3> 1 and d 7 <C 1) 
justifies theoretically why metastable states with a constant 
rate of increase of the q-entropy appear in the nonlinear dy- 
namical system. They argue that this law is reflected upon 
a transient behaviour in which the growth of the i.d.v. is 
almost a power law: w(t) cx t p ~< , with p 7 > 1 (and the 
exponent p 1 can be associated to a g-exponent via the rela- 
tionship p 7 = 1/(1 — q)). In fact, they introduce a method 
to compute the g-exponent and use it as a VIC (the APLE) 
in order to distinguish regular orbits from nearby weakly 
chaotic orbits. 

For an A-dimensional Hamiltonian H consider a par- 
titioning of the 2A-dimensional phase space S into a large 
number of volume elements of size 8 2N for some small S and 
let x(0) be the initial condition of an orbit located in a par- 
ticular volume element. The authors introduce the APLE: 

P'y . ELS follows: 



In 



Pi = 



( |w(t)| 2 \ 



where |w(t)| 2 = YlT=i \\ w k(t)\\ 2 an d Wfc(t) is one of the 
m deviation vectors of an orthogonal basis {wfc(t)} of the 
tangent space to S at the initial point x(0). Every Wfc(t) has 
length greater or equal to S and ti is a transient initial time 
of orbit evolution. 

Therefore, in case of regular motion, p 7 tends to 1, while 
for chaotic motion it grows exponentially fast. Moreover, the 
APLE shows oscillations around a constant value for a tran- 
sient time interval before the exponential growth for weakly 
chaotic orbits. The length of the time interval is directly 
connected with the hyper bolicity of the local phase space. 

Refer to Fig. 2 from ILukes-Gerakopoulos et all (|2008l ) 
to see the general behaviour of the APLE for regular and 
chaotic orbits. 



3.4 The FMFT 

Th e detailed description of the FMF lffl has been published 
on ISidlichovskv &: Nesvorn^l l| 19971 ). However, we should 
mention the conceptual aspects of this technique. Since the 
main advantage of the FMFT is that it improves the com- 
putation of the amplitudes and the frequencies given by the 



3 There is a downloadable version of the FMFT 
in the personal web page of David Nesvorny: 
http : //www. boulder . swri . edu/^davidn/ 



FMA jLaskarlll990l; lLaskar et al.|[l993: lLaskarlll993l , Il99d : 
iPapaphilippou fc Laskadll996l . Il998l : lLaskarl 120031 ). we de- 
cided to introduce only a short description of the latter. 

Let H be a non-integrable Hamiltonian which is also 
non-degenerate: 



det 



with I £ 6 C the action-like variables and v the fre- 
quency vector. The KAM theory states that for sufficient 
small values of the perturbation parameter e, there exists a 
Cantor set Q e of v that satisfies the Diophantine condition: 

llk-i/ll > 



where K e and m are constants and |k| is the norm of k. 
Furthermore, the system still possesses invariant tori with 
linear flow (KAM tori). The solutions lie on these invariant 
tori and are given in complex variables through their Fourier 
series: 



(in, u) 



where the coefficients a rn {i') depend smoothly on the fre- 
quencies i>i, i — 1, . . . , N. If we keep all the angle-like vari- 
ables fixed (0 G n", n" is the A-dimensional torus, fixed 
to = Oo), we obtain a frequency map in B defined as: 



Fe a : B -»■ fl e ; I 



p 2 (*- 1 (0o,I)) 



(6) 



where p2 is the projection on il e {p2(4>, v) = v). The FMA 
numerically gives through an iterative process the frequency 
map J- (i.e. amplitudes and frequencies) defined over the 
whole domain B, which coincides, within the numerical pre- 
cision adopted, with the _Fg of Eq. (JS| in the set of KAM 
tori. The frequency map T is obtained seeking for the 
quasiperiodic approximations of the solutions over a finite 
time-span through a finite series of terms: 



J(m k ,u) 



Once we obtain the quasiperiodic approximation from 
Eq. (0, we are able to establish a correspondence between 
the action-like variables (I G B) and the rotation num- 
bers (or frequency v ector v g Q e ). This correspondence is 
the so- called FMA (|Laskadll99oT ). Regular regions of the 
phase space allow a very accurate estimation of real rota- 
tion numbers in a given time interval (for instance, tx 1/T 4 
using the MFT with the Hanning window, cx 1/T 2 using 
the MFT without the Hanning window, cx 1/T using the 
usual FFT techniques, where [— T : T] is set as the whole 
time-span), and should not differ from the rotation numbers 
obtained in another time interval, considering a certain ac- 
curacy. Strictly speaking, the frequencies are only defined 
on KAM tori, the FMA numerically estimates over a fi- 
nite time-span a frequency vector for any initial condition, 
though. If the estimates of the rotation numbers given by 
the FMA vary (above the required accuracy) between dif- 
ferent time periods, it implies that the corresponding KAM 
tori are destroyed. Thus, the study of this frequency map 
T (i.e. the FMA) regarding its constancy in time provides 
important clues about the dynamics of the system. 
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Fo r a co mprenhensive discussion of the FMA refer to 
lLaskarl l|2003h . 

In order to calculate the frequencies with the FMFT, we 
need to compute the equati ons of motion in the ScTS model. 
We use the Taylor method jjorba fe Zoull2005h that proved 
to be very convenient for the task in D12. The precision in 
the computation of the coordinates is 10~ 15 . The numerical 
com putation of the VICs was done using the DOPRI8 rou- 
tine (|Prince fc D ormand 198l|). This routine achieved better 
results than Taylor when we required the simultaneous in- 
tegration of the coupled system of equations of motion and 
variational equations in the ScTS model (see also D12). The 
conservation of the energy with the DOPRI8 routine was 
~ lO" 13 , 1(T 14 . 

All the computations in the forthcoming investiga- 
tion were done using the following configuration. Hardware: 
CPU, 2 x Dual XEON 5450, Dual Core 3.00GHz; M.B., 
Intel S5000VSA; RAM, 4GB(4xlGB), Kingston DDR-2, 
667MHz, Dual Channel. Software: gfortran 4.2.3. 



4 THE FMFT AS A GLOBAL INDICATOR OF 
CHAOS 

In this section, we compare the performances of the FMFT 
and the SEflLCE as global indicators of chaos. 

4.1 Accuracy in the estimation of the frequencies 

In order to analyse the nature of each orbit we apply the 
FMFT to the following function in each degree of freedom 
(hereinafter, d.o.f.): 

M*k) = i>f(t h ) + ii>f{tk) (8) 

with ipj (j — x,y, z) a complex function, where ipf (the real 
part of tpj) corresponds to the Cartesian coordinates (x, y 
and z) and ipf (the conjugate part of ipj) corresponds to 
the Cartesian velocities (v x ,v y and v z ). The tk with fc = 
1, . . . , n are the trajectory samplings. Every set of points ipj 
generates a quasiperiodic approximation and we keep the 
frequency Vj(^ 0) with the maximum amplitude. We call 
this frequency the fundamental frequency of the system in 
the j d.o.f. 

Theoretically speaking, if the fundamental frequencies 
correspond to a regular orbit, they should remain constant 
in time. Numerically speaking, the variation exists and it is 
due to the precision of the computation. Thus, in order to 
determine such a precision for our particular experiment, we 
take samples of regular orbits, according to both convergent 
values of the L I and the MEGNO (< ln(t)/t or < 2.01, 
respectively, see lMafnone et ail |2011al )). 

We are aware that samples obtained from the use of 
VICs can only be statistically rich in regular orbits. That 
is, the actual values of the VICs are reached only at infinite 
time and as we have just truncated approximations of those 
values, some of the orbits might be missclassified. In order to 
be as certain as possible about the dominant type of orbits 
from the samples, we follow their behaviour for 10 3 T C (_E) 
(see Section |3.1[) . In a galaxy like the one represented by 
the ScTS model, the crossing time is 0.5 u.t. If we suppose 
that the Hubble time is of the order of 10 3 crossing times 
(a high value indeed), then we are integrating the orbits for 



~ 234 and ~ 14 Hubble times for the energy surfaces —0.1 
(T c (-0.1) = 1.17 x 10 2 u.t.) and -0.7, respectively. There- 
fore, 10 3 T C (E) seems to be a lapse of time large enough to 
provide a physical meaningful characterization of the mo- 
tion of the orbits and reliable convergent values of the VICs 
for most of the orbits in the sample. Finally, for all practi- 
cal purposes, most of the orbits of the samples are regular 
orbits. 

The samples are taken in the p XQ — p Z(l and the xo — zq 
start spaces for four different energy surfaces, namely: —0.1, 
-0.3 (T c (-0.3) = 24 u.t.), -0.5 (T c (-0.5) = 12 u.t.) and 
—0.7 (i.e. a total of eight samples). There are 100 i.e. in 
each sample and we integrate the equations of motion for 
4.5 x 10 2 T C [E) for every i.e. Although 10 3 T C (E) is the time- 
span we use to compute the VICs reliably, the results in 
the determination of the fundamental frequencies with the 
FMFT are not very different using only 3 x 10 2 T C (_E). Be- 
sides, the CPU time is certainly reduced. Therefore, we ap- 
ply the FMFT on two 50% overlapping tim e intervals of 
3 x 10 2 T C {E) ()Wachlin fc Ferraz-Mellol Tl998l fl. The num- 
ber of points used in the estimation of the frequencies is 
8192 which makes ~ 27 points for each T C (E). Finally, we 
estimate the fundamental frequencies v^ 1 for each time in- 
terval and calculate the differences (the accuracy in the de- 
termination of the frequencies sh ould be ~ 10 -9 according 
to Section T3.4I and [Muzziol (|2006l )). The concomitant results 
for the energy surfaces —0.1 and —0.7 for both the p XQ — p ZQ 
and the xq — zq start spaces are shown in Fig. [2] 

In Fig. [2] we show that the differences in the computed 
values of the fundamental frequencies between both over- 
lapping time intervals for regular orbits (actually, some of 
them might be chaotic) is not generally below the theoretical 
value, particularly in the p XQ — Pz start space. Thus, there 
exists a range of values where the regular orbits could not 
be reliably distinguished from the chaotic orbits by means 
of the computation of the fundamental frequencies with the 
FMFT. 

In Table [3] we show the corresponding arithmetic means 
and standard deviations of the variations of the fundamen- 
tal frequencies for the energy surfaces: —0.1, —0.3, —0.5 and 
—0.7 on the p XQ — p Z() start space. The arithmetic means 
show an average conservation between the 7th and the 8th 
decimal in the fundamental frequencies. Nevertheless, the 
standard deviations oscillate around an average of the 5th 
decimal. Therefore, there is an interval where the FMFT can 
not distinguish clearly between regular and chaotic orbits in 
the experiment. However, we can infer certain limits which 
can help us to identify if the orbit is regular or chaotic us- 
ing the FMFT. For instance, the variations in the computed 
fundamental frequencies for the regular orbits only affect 
decimals higher than the 4th (see Fig. Thus, if there 
is a difference in the computation of the fundamental fre- 
quencies of the orbit in decimals lower than the 4th, it is 
highly probable that the orbit is chaotic. This resu lt is in 
complete agreement with the one shown in iMuzziol (2006) 
for the same ScTS model an d with the practical l imit used 
also in lAquilano et all j2007l ): lMuzzio et all (|20Qgh . 



4 The time intervals are [0 : 3 X 10 2 T C {E)\ and [1.5 X 10 2 T C (E) : 
4.5 x W 2 T C (E)]. 
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Figure 2. Variations of the fundamental frequencies (in logarithmic scale) computed with the FMFT for four samples of 100 regular 
orbits on each d.o.f.: on the p xg — p ZQ start space and on the energy surfaces —0.1 (top left panel) and —0.7 (top right panel), and on 
the xo — zo start space and on the energy surfaces —0.1 (bottom left panel) and —0.7 (bottom right panel). The missing points are due 
to the computational precision. 



Table 3. Arithmetic means and standard deviations of the differ- 
ences on the fundamental frequencies between both overlapping 
time intervals for the samples of regular orbits used in Fig. [2] The 
values are ordered by energy surface and d.o.f. 



Energy —0.1 


X 




V 


2 




Arithmetic means 


~ 1.2 - 


9 


2.3~ 7 


~ -3.4 


-7 


Standard deviations 


~ 1.2- 


7 


~ 2.5~ 5 


- 1.7- 


5 


Energy —0.3 


Arithmetic means 


~ -1.8 


-8 


~ 1.2" 6 


- 5.9- 


(i 


Standard deviations 


~ 4.6" 


6 


~4.3~ 4 


- 5.9- 


4 


Energy —0.5 


Arithmetic means 


~ 5.3" 


9 


~ -1.9- 8 


~ 1.1- 


7 


Standard deviations 


~ 5.3" 


7 


~ 2.7" B 


~ 1.3- 


5 


Energy —0.7 


Arithmetic means 


~ 3.5" 


7 


~ -6.9" 8 


~ -1.9 


-S 


Standard deviations 


~ 1.6" 


5 


~1.3- 5 


- 1.2- 


3 



The results for the xo — zq start space are similar than 
those for the p XQ — p ZQ start space. Nevertheless, for lower 
energy values (Fig. bottom right panel), the averaged pre- 
cision in the computation of the frequencies improved sub- 
stantially. 



We also make some experiments with the FMFT, using 
synthetic signals for verifying the precision in the compu- 
tation of the fundamental frequencies. The results were re- 
markable bett er and are closer to the theoretical estimates 
<|Muzzidl20(J6h . Nevertheless, this is not the case for the rel- 
atively complex ScTS model and the precision in the esti- 
mations is clearly diminished. Thus, we use the same model 
and the same spaces we study to calibrate the FMFT. 



4.2 A comparative evaluation between the FMFT 
and the SE/LCE 

We are interested in testing the performance of the FMFT as 
a global indicator of chaos. Thus, we apply it together with 
a variational indicator in order to characterise the phase 
space portraits of two different regions. One of the regions 
is located in the p XQ ~Pz start space (the region is identified 
in the top panel of Fig. [3]), whereas the other is located in 
the xo — zo start space (which is identified in the bottom 
panel of the same figure). Both regions belong to the energy 
surface -0.7. The VIC selected for this study is the SE/LCE 
(see Section f3.ip . 

In order to apply the FMFT as a global indicator of 
chaos, we follow the same procedure as described in Section 



14.11 Thus, we consider two overlapping time intervals of 10 
periods. The first time interval is [0 : 7 x 10 3 ] u.t. and the 
second time interval is [3.5 x 10 3 : 10.5 x 10 3 ] u.t. According 



to Section T3. 41 the accuracy in the determination of the fre- 
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Figure 3. On the top panel we show in a red square the region we 
will study in the p XQ — p ZQ start space and on the bottom panel, 
the one we will study in the xq — zg start space. Both regions are 
located in the energy surface defined by the constant value —0.7. 
We take a total integration time of 7 X 10 3 u.t. The phase space 
portraits are generated in grey-scale with the MEGNO. Values 
of the MEGNO below the threshold 2.01, i.e., regular orbits, are 
depicted in white. In grey we enclosed the values of the MEGNO 
in the interval (2.01 : 20), which means moderate to strong chaos. 
In black, we show values of the MEGNO above 20, i.e. regions 
with very strong chaos. 

quencies should be ~ 10 -11 . The computation is over 624100 
i.e., within the p XQ — p ZQ start space, covering the region de- 
picted on the top panel of Fig. [3] and over 596258 i.e., within 
the xo — zo start space, covering the region depicted on the 
bottom panel of the same figure. 

To use the FMFT as a global indicator of chaos and 
compute the correspon ding phase space portraits, we n eed to 
define a new quantity (jWachlin fc Ferraz-Me llo 1998). Such 
a quantity is log(AF), AF being |i4 1} - *4 2) | + \vy ] - vf> | + 
, where is the frequency (computed with the 
FMFT) associated with the d.o.f "j" for the time interval 
"(i)". All the fundamental frequencies have to be computed 
in both intervals. Finally, the phase space portraits shown by 
the quantity log(Ai ? ) consist of 622521 i.e. for the p XQ — p ZQ 
start space and 594690 i.e. for the xq — zq start space. 

In order to make appropriate comparisons on each start 
space, we need to evaluate the log(Af) and the SE/LCE on 
equal standings. Thus, it is necessary to decide a cut-off 
level for the distributions of their values to use the same 
scales for both indicators. 

In Fig. [4] we show the histograms according to the final 
values of the techniques. The histograms on the top panel 



1e+06 




-10 -S -S A -2 

Global indicator of chaos 




Global indicator of chaos 

Figure 4. Histograms in logarithmic scale showing the bimodal 
distribution of both global indicators of chaos on the p xg — p ZQ 
(top panel) and the xq — zq (bottom panel) start spaces on the 
energy surface —0.7. The cut-off level of 10 -8 is depicted on the 
top panel. 



correspond to the p XQ — p zg start space and the histograms 
on the bottom panel correspond to the xo — zo start space. A 
clear bimodal distribution (one peak corresponding to reg- 
ular orbits, and the other to chaotic orbits) is shown for 
both indicators. In case of the p X(l — p ZQ start space, the in- 
dicators have less than 1% values for the orbits below 10~ 8 
(pointed out on the top panel of the same figure). Thus, the 
scales on the p XQ — p za start space go from 10 -8 to 0. The 
similarities between both histograms on the p XQ — p ZQ start 
space make easier the selection of a common cut-off level. 
However, there are important differences between the his- 
tograms on the xo — zo start space. The SEiLCE shows a 
bimodal distribution within the interval [~ 10~ 10 : 0) with 
its peaks well separated. The peaks in the distribution on 
the xo — zo start space match the ones found for the same 
indicator, but on the p X(j ~p ZQ start space (see both panels in 
Fig. [J). This is not the case for the log(AF). The peak that 
corresponds to the regular component is beyond the value 
10~ 10 (see Fig. [4j bottom panel and Fig. [2] bottom right 
panel). Nevertheless, we find that the number of orbits with 
values of the SE/LCE below the value 10 -10 are negligible 
and the log(Ai ? ) does not show further structure below that 
value. Thus, we decided to use that value as the cut-off level 
in the distributions on the xo — zo start space to build the 
scales. 

Notice also that, according to log(AF) on the p XQ — 
p ZQ start space, the peak corresponding to regular orbits is 



© - RAS, MNRAS 000,[THT9] 



10 N. P. Maffione et al. 





0.7 
0.6 
0,5 
0.4 
0.3 
0.2 
0.1 






-1 

-2 
-3 
-4 
-5 
-6 
-7 
-B 
-9 
-10 




0.4 0.6 



Figure 5. Phase space portraits corresponding to both regions within the p xg —Pz (top panels) and the xq — zq start (bottom panels) 
spaces on the energy surface —0.7. The samples used to compute the SEZLCE consist of 624100 and 596258 (top and bottom left panels, 
respectively) i.e. and 7 X 10 3 u.t. The samples used to compute the log(AF) consist of 622521 and 594690 (top and bottom right panels, 
respectively) i.e. and two overlapping time intervals of 7 X 10 3 u.t. each. Chaotic regions are identified by the values of the indicators 
that are higher than ~ —4, ~ —3 (notice that the values of the indicators are in logarithmic scale). See text for further details. 



~ 3.1 x 10 . This value confirms the conclusions drawn 
in Section 14.11 with a sample composed only of 100 orbits 
for the same start space. That is, if we use the arithmetic 
means for the energy surface —0.7 (Table [3]), we arrive to 
a similar value of ~ 1.9 x 10 . Furthermore, the minimum 
between both peaks (which gives a clear distinction between 
the regular and the chaotic component) in the p X(l — p ZQ 
start space is taken around the value 10 -4 . This value also 
matches the one considered to separate regular from chaotic 
orbits for the same start space in Section H. II 

In Fig. [5]we present the phase space portraits given by 
the SEZLCE (left panels) and the log(A.F') (right panels) for 
both selected regions in the p XQ — p ZQ (top panels) and the 
xo — 2o (bottom panels) start spaces. 

The phase space portraits of the region inside the 
Px —Pz start space given by the SEZLCE and the \og(AF) 
present the same structures. In spite of this high level of co- 
incidence in the description of the p XQ —p ZQ start space, there 
are differences in the xq — zo start space. On the bottom right 
panel of Fig. [5] the quantity log(AJ ? ) shows an amount of 
spurious structure^] not shown by the variational indicator 
(bottom left panel of the same figure). In fact, the spurious 
structures make difficult to select a threshold to distinguish 

5 Some of them, are due to the phenomenon of Moire, which 
is inher ent in discrete Fourier transform techniques, Barri o et al.l 
(20093). 



between regular and chaotic orbits. Still, it is possible to 
define a cut-off value for the log(A.F), but the procedure is 
more handmade. For example, if we take as regular orbits 
those preserving > 4 decimals in their fundamental frequen- 
cies (see Section al I and the previous discussion), we recover 
the portrait given by the SE/LCE for the xo — zo start space. 
However, the SE/LCE (as well as most of the VICs studied 
in previous papers, for instance Mil and D12) shows a large 
separation of the different kinds of motion in both portraits 
(see also Fig. [4]), and the choice of a threshold (if it is not 
already defined) is easier. 

These results, which are not as clear as those given 
by the variational indicators in terms of motion separation, 
make the FMFT a second choice to study the global dynam- 
ics in a divided phase space. 

The procedure used to determine the chaoticity or reg- 
ularity of the orbits with the FMFT is standard. The defi- 
ciency in the descriptions is basically due to the method's 
sensitivity to its parameters. Therefore, the reliability of the 
FMFT as a global indicator of chaos is, comparatively speak- 
ing, limited. 

4-2.1 The computing times 

The VICs need to compute the equations of motion together 
with the variational equations to determine the nature of the 
orbits. This process may take a large amount of time. On 
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the other hand, the FMFT works with the computation of 
the frequencies, which is a fast process, but the determina- 
tion of the log(AF) to characterise the nature of the orbits 
involves further calculations. For example, the computation 
of the SE/LCE (a fast VIC) to obtain the top and bottom 
left panels of Fig. retook ~670 hs and ~330 hs, respectively, 
by using the DOPRI8 routine. Instead, the computation of 
the fundamental frequencies (nothing more than the the line 
of biggest amplitude in each d.o.f. was computed) with the 
FMFT for 624100 and 594690 i.e. took ~30 hs for both ex- 
periments, after the integration of the equations of motion. 
The time spent in this integration for 7 x 10 3 u.t. was of 
~570 hs and ~280 hs for the regions in the p XQ — p Z(l and 
the xo — zo start spaces, respectively. In both cases we used 
the Taylor method (see Section l3~4")) . Furthermore, if we want 
to use the FMFT as a global indicator of chaos, we need to 
compute the log(AF), i.e. we need to compute the frequency 
vector for two time intervals. If we overlap the time intervals 
in a 50%, we economize on computing time. The total time 
interval must be 1.05 x 10 4 u.t. under these circumstances. 
That is, 50% larger than the one used by the variational 
indicators. Finally, the time required by the FMFT extends 
from ~600 (~570 + ~30 hs) to ~885 hs and from ~310 
(~280 + ~30 hs) to ~450 hs to obtain the portraits on the 
top and bottom right panels of Fig. [5] respectively. 

The advantage of the fast computation of the frequen- 
cies is certainly lost in the time-consuming process involved 
in the computation of the log(AF'). 

In conclusion, the FMFT as a global indicator of chaos 
seems to be rather inconvenient when there are other alter- 
natives such as fast VICs like the SE/LCE. 



5 THE SELECTION OF PROPER VICS 

The aim of this section is to select the proper VICs to be 
used in the forthcoming experiments on the p Xo — p ZQ and 
the xo — zo start spaces of the ScTS model (Section [6]). 

We first study a small region in the p XQ — p ZQ start 
space (within the intervals p XQ £ [0.87 : 1.02] and p Zo £ 
[1.39 : 1.54]) on the energy surface —0.7, shown in the red 
square in Fig. [6]). We use the MEGNO, the SE/LCE, the 
FLI, the OFLI and the APLE. This first experiment pro- 
vides us with a comparison between quite new techniques 
(the SE/LCE and the APLE) and VICs with excellent per- 
formances shown in previous studies (the MEGNO and the 
FLI/OFLI, e.g. in Mil and D12). 

In addition, we use the same 10 3 characteristic times 
along the experiments. 

We are interested in the resolving power and the speed 
of convergence of the VICs. Therefore, we compute 10201 
i.e. within the region marked in Fig. [5] for 7 x 10 3 u.t. 

5.1 The resolving power of the VICs 

First, we briefly discuss the resolving power of the VICs 
under the circumstances of the experiment. 

As we are looking for VICs to study big samples of i.e., 
we consider the final values of the indicators rather than 
their time evolution curves. Furthermore, when the indicator 
grows or decreases exponentially, it is convenient to stop the 
integration process at a particular saturation value. Thus, 




0.5 1 1.5 2 2.5 



Figure 6. We signal with a red square the region in the p xg — p Z(t 
start space (energy: —0.7) under analysis. The grey-scale portrait 
is obtained with the final values of the MEGNO. In white, we 
identify the regular orbits and thus the MEGNO values < 2.01. 
In grey, we identify regions of moderate to strong chaos, that is 
MEGNO values in the interval (2.01 : 20). In black, with MEGNO 
values > 20, we identify very strong chaos. 



instead of simply registering the final values of the indica- 
tor, we can record the time nee ded to reach such sa turation 
value, i.e. the saturation times l|Skokos et al.ll2007l ). The fi- 
nal values together with the saturation times of the VIC 
should be considered for the description of the phase space 
portrait (Mil and D12). 

In Fig. [7] we present the phase space portraits given by 
the final values of the MEGNO and the SE/LCE (top left 
and top right panels, respectively); the OFLI final values 
(middle left panel) and the saturation times (middle right 
panel); and the APLE final values (bottom left panel) and 
the saturation times (bottom right panel). The phase space 
portraits given by the FLI are very similar to those shown 
with the OFLI and they are not included for the sake of 
simplicity. 

The region under analysis has a regular component (the 
main central structure) and a chaotic component (surround- 
ing the main central structure). Moreover, there are two dif- 
ferent regions of chaos. The region that surrounds the reg- 
ular component is composed of very strong chaotic orbits 
while the band on the right side of every panel in Fig. [7] is 
filled with moderate to strong chaotic orbits. 

The values of the MEGNO and the SE/LCE do not 
grow exponentially for chaotic orbits as it happens, for in- 
stance, with the OFLI (and the FLI) or the APLE. Then, a 
saturation value is not defined and a saturation time is not 
recorded in the computing process. However, the values of 
the MEGNO and the SE/LCE can be very large for strong 
chaotic orbits and very different from their values for regu- 
lar and mild chaotic orbits. In a divided phase space with 
regions of very strong chaos, the large values of these indica- 
tors for strong chaotic orbits may hide the small differences 
between regular and mild chaotic orbits. Indeed, the differ- 
ences among strong chaotic orbits are generally not as im- 
portant as the differences between strong and mild chaotic 
orbits or between mild chaotic orbits and regular orbits. 
Thus, a saturation value for the MEGNO or the SE/LCE 
is also very useful in this kind of scenario. As shown in Fig. 
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Figure 7. We show the phase space portraits of the different VICs using a total integration time of 7 X 10 3 u.t. for the region displayed 
on Fig. [6] We present the final values of the MEGNO (top left panel); the final values of the SEZLCE (top right panel); the final values 
(middle left panel, in logarithmic scale) and the saturation times (middle right panel) of the OFLI; and the final values (bottom left 
panel) and the saturation times (bottom right panel) of the APLE. The cold colours correspond to regular regions while the warm colours 
correspond to chaotic domains for the final values. The colour scale is inverted for the saturation times. 



[6] we take the saturation value 20 for the MEGNO and its 
associated value ~ 0.0057 for the SEZLClfi 

In case of the OFLI (and the FLI) or the APLE, the 
chaotic orbits grow exponentially fast and a saturation value 
is needed to avoid the propagation of errors in their compu- 
tation. The saturation value in Fig. [7] for the OFLI (middle 
panels) is ~ ln(10 20 ) (see Section fjOJ- The saturation value 
in Fig. for the APLE (bottom panels) is ~ ln(10 20 )/ ln(£) 
(see Section 133)) . 

6 We can estimate the saturation value (x-y) for the SEZLCE from 
the expression: 20 ~ X7/2 X t, with t = 7 X 10 3 u.t. (see Section 

ED - 



In the top panels of Fig. [7J we can see the regular and 
the chaotic components well separated by the MEGNO (left 
panel) and the SEZLCE (right panel). Furthermore, the two 
different chaotic regions are clearly shown: the region of 
strong chaos (surrounding the central structure with regular 
orbits) and the band on the right side of the panels filled with 
moderate to strong chaotic orbits. There is no structure at 
all in the main central structure composed of regular orbits 
(top left panel) because of the MEGNO 's asymptotic thresh- 
old (Section 13. ip . Also, there is no further structure in the 
strong chaotic region due to the selected saturation value. 
The SEZLCE shows a similar portrait (top right panel). 

In middle panels of Fig. [7] we use the OFLI to describe 
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the region under analysis. On the one hand, we recover 
the phase space portraits (middle left panel) given by the 
MEGNO and the SEZLCE in the top panels of the same 
figure. On the other hand, the band on the right side of 
both panels shows a more detailed structure with the OFLI. 
This is due to the fact that the saturation value taken for 
the OFLI covers a wider range of chaotic orbits than the 
values chosen for the MEGNO or the SEZLCE. The satu- 
ration times (middle right panel) enhance the information 
given with the final values of the OFLI, specially within the 
strong chaotic region. 

In the bottom panels of Fig. we present the phase 
space portraits given by the APLE final values (left panel) 
and saturation times (right panel). The description is very 
similar to those given by the other VICs if we consider both 
the final values and the saturation times of this indicator. 
For example, the identification of the band on the right side 
of the panels needs the information of the saturation times. 
The saturation times show a clear distinction between the 
region of strong chaos surrounding the main central struc- 
ture and the region with mild and strong chaotic orbits in 
the band on the right side of the panels. This is not the case 
for the other VICs (the MEGNO, the SEZLCE, the OFLI 
and the FLI) for which the saturation times play a com- 
plementary role and the final values were enough for such 
identification. Nevertheless, the degree of coincidence in the 
descriptions of the four VICs (plus the FLI) seems to be 
high. Thus, there is no decisive advantage in favour of one 
of them. 

5.2 The speed of convergence of the VICs: 
thresholds 

Now we focus on another important characteristic of the 
VICs which is fundamental for an efficient study of large 
samples of i.e.: the speed of convergence. 

We use the same total integration time (7 x 10 3 u.t.) 
and the same i.e. (10201) used in Section 15.11 The initial 
10 2 u.t. of the time interval does not give reliable informa- 
tion for our purposes and it is not considered in this ex- 
periment. The remaining 6.9 x 10 3 u.t. are divided into 70 
sub-intervals of ~ 98.57 u.t. each. We compute the num- 
ber of chaotic orbits in each sub-interval. Before that, we 
had to identify the chaotic orbits. So, we needed accurate 
thresholds for each VIC. Many of them have theoretical or 
empirical estimations of their thresholds. Yet, they are only 
approximations and should be adjusted to each particular 
situation. Hence, we generalize the experiments shown in 
Mil and D12 and consider a wide range of thresholds for 
the different VICs. Finally, for each sub-interval, we com- 
pute the number of chaotic orbits for every VIC according 
to the different thresholds. 

Performing several tests (i.e. considering different 
thresholds and VICs), we find that a percentage of ~ 83% 
for the chaotic component is stable and prevails in the ex- 
periment. We take it as the "true" percentage for the chaotic 
component. Having such a "true" percentage, we can adjust 
the thresholds of the VICs under consideration. The fastest 
convergency of the VIC to a stable percentage of chaotic 
orbits close to the "true" percentage gives the most appro- 
priate threshold for that particular indicator. 

In order to identify such appropriate thresholds for ev- 



ery VIC under consideration, we define the "selected thresh- 
olds". These thresholds allow the associated indicator to 
reach the "true" percentage of chaotic orbits within the total 
integration timqj. 

In Fig. [8] we present the variation of the selected thresh- 
olds with the integration time for the APLE (top left panel), 
the MEGNO (top right panel), the FLI (bottom left panel) 
and the OFLI (bottom right panel). The SEZLCE is not in- 
cluded in the experiment because it has not a defined mech- 
anism to compute its threshold. The filled circles in Fig. [H] 
indicate the selected thresholds that allow the convergency 
of the corresponding VIC to a stable, close to the "true" 
percentage of chaotic orbits. The filled squares in the same 
figure indicate the selected thresholds that offer the fastest 
convergency. In other words, the filled squares indicate the 
most appropriate thresholds for a given indicator. 

The APLE has the theoretical time- 
independent threshold ~ 1 (see Section 13.31 and 
iLukes-Gerakopoulos et all |2008) for further details). 
However, we consider threshold values between ~ 1.2 and 
~ 1.5. The threshold values below ~ 1.2 characterise 100% 
of the orbits as chaotic orbits for the whole time interval 
and thus, they are discarded. The threshold values above 
~ 1.5 do not arrive at the "true" percentage of the chaotic 
component within the time interval and thus, they are also 
discarded. If we use the threshold values in the interval 
~ [1.4 : 1.5], the APLE converges to stable percentages 
of the chaotic component close to the "true" percentage 
(see top left panel of Fig. |HJ). If we use the threshold 
values around 1.4, the APLE's convergency is the fastest 
(approximately by 5.5 x 10 3 u.t.). 

The MEGNO, like the APLE, has also a time- 
independent thresh old: ~ 2 (see Section 13.11 and 
ICincotta et alj (|2003T ) for further details) . For the same rea- 
sons given for the APLE, we consider threshold values above 
the theoretical approximation, i.e. between ~ 2.8 and ~ 4. 
The threshold values below ~ 2.8 characterise 100% of the 
orbits as chaotic orbits for the whole time interval and thus, 
they are discarded. The threshold values above ~ 4 are much 
larger than the theoretical asymptotic estimation of ~ 2 
and thus, they are also discarded. If we use the threshold 
values close to ~ 4, the MEGNO reaches percentages of 
the chaotic component close to the "true" percentage more 
rapidly (see top right panel of Fig. |SJ). Nevertheless, using 
this range of selected thresholds the percentages are not con- 
vergent. Moreover, the percentages oscillate within the in- 
terval ~ [84% : 87%] depending on the threshold selected. 
In other words, the MEGNO gives a higher percentage of 
chaotic orbits than the prevailing percentage shown by the 
other VICs. This experiment confirms the results yielded in 
Mil and D12. 

The FLI and the OFLI have both similar behaviour and 
time-dependent thresholds to distinguish chaotic motion 
from regular motion (see Section 13.21 and iFouchard et al.l 
(2002) for further details). After computing the logarithm 
in their definitions, a tentative threshold for them is In(t) 
with t the integration time. We consider the general expres- 



7 Actually, the selected thresholds allow the associated VIC to 
reach a percentage of chaotic orbits between 82-84% within the 
total integration time. 
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Figure 8. We present the selected thresholds for the APLE (top left panel), the MEGNO (top right panel), the FLI (bottom left panel) 
and the OFLI (bottom right panel). The sample consists of 10201 i.e. for 7 X 10 3 u.t. in the region shown in Fig. [6] See text for further 
details. 



sion a x \n(t) in order to study a variation of the thresholds 
changing the values of the parameter a. 

The values of the parameter a for the FLI belong to the 
interval ~ [1.2 : 1.6]. The theoretical threshold with a ~ 1 is 
not efficient because it does not arrive at stable percentages 
of the chaotic component close to the "true" percentage for 
the whole time interval. The same applies to the thresholds 
with a < 1.2 so they are discarded. Furthermore, the thresh- 
old values with a > 1.6 do not reach at the "true" percentage 
of the chaotic component within the time interval and thus, 
they are also discarded. If we use the threshold values in the 
interval ~ [1.2 : 1.6] in the time interval ~ [5.7xl0 3 : 7xl0 3 ] 
u.t., the FLI converges to stable percentages of the chaotic 
component close to the "true" percentage (see the bottom 
left panel of Fig. [SJ). If we use the threshold values with 
a ~ 1.3, the FLI's convergency is the fastest (approximately 
by 5.7 x 10 3 u.t., i.e. 1.3 x ln(5.7 x 10 3 ) ~ 11.2 in the figure). 
These results are very similar to those obtained earlier with 
the APLE. Yet, there is an important difference: the thresh- 
old selected for the FLI (a ~ 1.3) is closer to the theoretical 
approximation (a ~ 1) than the threshold selected for the 
APLE (a ~ 1.4). Therefore, this indicates that the FLI may 
have a better theoretical approximation of its threshold than 
the APLE. 

Finally, we also observe similar results for the OFLI. 
However, the parameter a varies from ~ 0.5 to ~ 1.5, i.e., 
for the first time the interval of selected thresholds includes 
the theoretical estimation (a ~ 1). That is, if we use the 
theoretical approximation, the OFLI reaches a stable per- 



centage of chaotic orbits close to the "true" percentage of 
83%. Nevertheless, the convergency arrives at the end of 
the time interval, i.e. ~ 6.9 x 10 3 u.t. The threshold val- 
ues below ~ 0.5 characterise 100% of the orbits as chaotic 
orbits for the whole time interval and thus, they are dis- 
carded. The threshold values above ~ 1.5 do not arrive at 
the "true" percentage of the chaotic component within the 
time interval and thus, they are also discarded. If we use the 
thresholds with a in the interval ~ [1 : 1.4] and within the 
time interval ~ [5.7 x 10 3 : 7 x 10 3 ] u.t., the OFLI converges 
to stable percentages of the chaotic component close to the 
"true" percentage (see the bottom right panel of Fig. [8}. If 
we use the threshold values with a ~ 1.07, the OFLI's con- 
vergency is the fastest (approximately by 5.7 x 10 u.t., i.e. 
1.07 x ln(5.7 x 10 3 ) ~ 9.2 in the figure). These results are 
very similar to those obtained with the FLI as mentioned 
before. Yet, once again there is an important difference: the 
threshold selected for the OFLI (a ~ 1.07) is closer to the 
theoretical approximation (a ~ 1) than the threshold se- 
lected for the FLI (a ~ 1.3). Therefore, this indicates that 
the theoretical approximation of the thresholds of the OFLI 
may work better. 



© - RAS, MNRAS 000,[THT9] 



Methods in chaos detection 15 



6 GLOBAL ANALYSIS OF THE START 
SPACES 

In this section we proceed to study some dynamical aspects 
of the ScTS model by means of the complementary use of 
the VICs selected in Section [5] 

We consider samples of around one millon i.e. on the 
Px — Pz and the xo — 20 start spaces for four differ- 
ent energy surfaces. Even though we do not study the re - 
maini ng five start spaces indicated in ISchwarzschiidl (|l993h : 
iPapaphilippou k, Laskarl (|l998l ') . we gain sufficient informa- 
tion in order to briefly discuss some main features of the 
dynamics of the system. 

The total integration times extend to 1.17 x 10 5 u.t. for 
the energy surface —0.1 (Section 14. 1|) . The CPU times be- 
come critical and will be considered a priority in the selection 
of the VICs. The FLI family of indicators presents the most 
versatile indicators, showing good performances in all the ex- 
periments (see also Mil and D12 in which they are included 
in the corresponding CIsF). The MEGNO is also another 
candidate because of its good performance (also included 
in the CIsF of Mil and D12). The SEZLCE seems to be a 
good alternative to the MEGNO when studying big samples 
of orbits. Furthermore, the MEGNO and the SEZLCE show 
a lower computational cost in the experiment than that of 
the FLI family of indicator^]. Therefore, we select the LI 
(the least time-consuming indicator given a fixed total inte- 
gration time, see e.g. D12) and the pair MEGNO/SEZLCE 
as the VICs to study the ScTS model. 

6.1 The p xg — Pz start space 

We have already mentioned that the SEZLCE has certain ad- 
vantages when we study big samples of orbits (see Section 
13. 1|) . Nevertheless, it does not have any defined procedure 
to determine a threshold in order to distinguish between 
regular and chaotic motion. Indeed, this fact is a serious 
problem if we want to separate regular motion from weak 
chaotic motion. We need a threshold for the SEZLCE to cor- 
rectly describe the p xa — p ZQ start space of the ScTS model 
through the final values of the indicator. For the determina- 
tion of such a threshold we require again a "true" percentage 
of chaotic orbits for the sample under study. The MEGNO 
shows high sensitivity to changes in its threshold (Section 
15.21 Mil and D12) and hence it is not a reliable indicator 
to calibrate other VICs. Besides, the LI has a theoretical 
threshold (ln(t)/t, with t the total integration time), which 
is a still valid estimation as a decent first approximation. We 
start with the corresponding theoretical estimation^] and 
calibrate them using inspections of different profiles of i.e. 
within the sample. Then, we find reliable thresholds for the 
LI for every energy surface. We compute the "true" percent- 
ages of chaotic orbits with the LI and calibrate the thresh- 
olds of the SEZLCE to best fit these values. In Tabled we 

8 In this tes t on computing times, we consider the definition of the 
FLI given in lFroeschle &: Legal 1 120001 ) , where the authors compute 
the evolution of the FLI using one deviation vector. 

9 The theoretical estimations: 9.9743 X 10 -5 (u.t.) -1 for the en- 
ergy surface -0.1; 4.20242X10" 4 (u.t.)" 1 for -0.3; 7.82722xl0~ 4 
(u.t.)" 1 for -0.5 and 1.264809 X 10~ 3 (u.t.)" 1 for -0.7 give over- 
estimated percentages of the chaotic component. 



present the results ordered by energy surface. Notice that 
the thresholds associated with the SEZLCE are considerably 
smaller than the thresholds of the LI, because the estima- 
tion of the 1LCE given by the slope of the MEGNO converges 
more rapidly to the 1LCE than the st andard procedure (see 
Section [3~T1 and ICincotta et all (|2003h V 

A comment regarding the efficient way to work with the 
LI and the SEZLCE is in order before going forward. On one 
hand the LI is a reliable, but slow indicator so that it should 
be applied for longer integration times but on smaller sam- 
ples in order the computing times to be moderate. Then, 
we can estimate the percentages of chaotic orbits. On the 
other hand, we have the SEZLCE, which has a fast speed 
of convergence so it can be applied for shorter integration 
times but on the whole sample of orbits (once again, the 
computing times are reasonably short). Finally, we calibrate 
the thresholds of the SEZLCE adjusting the percentages of 
chaotic orbits to those obtained with the LI. This is an ef- 
ficient way to take advantage of the different characteristics 
(reliability and speed of convergence) of both VICs. Never- 
theless, the priority in this experiment is to obtain reliable 
percentages of chaotic orbits rather than saving computing 
time in order to obtain the more accurate phase space por- 
traits. This is the reason why we previously applied both the 
LI and the SEZLCE on the same samples and for the same 
total integration times. 

Now that we have the thresholds for the SEZLCE, we 
apply the latter to a sample of 1000444 i.e. to study the 
Px —Pz start space for the following four energy surfacej 1 ^!: 
-0.1, -0.3, -0.5 and -0.7 (see Fig.©. 

In Fig. [9] we observe that independently of the energy 
surface, the chaotic component is the dominant one on the 
phase space portraits (yellow and orange colours. See also 
column three and five in Table |4| . The regular component 
(red, blue and black colours) increases for lower energy val- 
ues. However, this increment is not considerable yet. 

On the top left panel of Fig. [5] we observe a clear separa- 
tion of both components. We show a fully connected chaotic 
component for values of p XQ < 1.7 and the regular compo- 
nent for values of p XQ > 1.7, besides some hyperbolic struc- 
tures. 

These structures grow as we move to lower energy sur- 
faces. They are resonances that start to overlap with each 
other and occupy the regular component. We can also dis- 
tinguish a separation inside the chaotic component. The 
fully connected chaotic domain moves back to lower val- 
ues of p XQ . The region is occupied by another chaotic do- 
main characterised by a resonance overlapping regime and 
a lower Lyapunov number (top right panel of Fig. On 
the energy surface —0.5 (bottom left panel of Fig. ©, the 
fully connected chaotic domain has noticeably moved back- 
wards and the resonances occupy extended regions of the 
Px a — Pz start space. Furthermore, a major resonance is 
seen around the value p X() ~ 0.5 and minor resonances are 
seen around the values p XQ ~ 1 and p Xf> ~ 1.3. On the bot- 



10 We use a cut— off level around 10" 9 to have the same scales in 
Fig. [9] for every energy surface. Then, the percentages of orbits 
that we do not consider are negligible. Those percentages are 
~ 0.77% in the energy surface -0.1, ~ 0.06% in -0.3, ~ 0.03% 
in -0.5 and ~ 0.02% in -0.7. 
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Table 4. For every energy surface on the p xo —pz start space, we present the threshold associated to the LI and the estimated percentage 
of chaotic orbits; the threshold associated to the SEZLCE and the estimated percentage of chaotic orbits. 



Energy surface 


Threshold (LI) 


Chaotic orbits (%) -LI 


Threshold (SEZLCE) 


Chaotic orbits (%) -SEZLCE 


-0.1 


1.2 x 10~ 4 


~ 71.46% 


1.4 x 10" 5 


~ 66.83% 


-0.3 


5.8 x 10~ 4 


~ 65.89% 


6.4 x 10" 5 


~ 65.99% 


-0.5 


1.25 x 10~ 3 


- 61.79% 


1.3 x 10~ 4 


~ 63.05% 


-0.7 


1.7 x 10" 3 


~ 65.09% 


2.7 x 10~ 4 


~ 62.87% 
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Figure 9. Phase space portraits of the SEZLCE for a sample of 1000444 i.e. on the p XQ —pz start space and for four energy surfaces. Top 
left panel: phase space portrait corresponding to the energy surface —0.1 and using a total integration time of 1.17 X 10 5 u.t. Top right 
panel: —0.3 and 2.4 X 10 4 u.t. Bottom left panel: —0.5 and 1.2 X 10 4 u.t. Bottom right panel: —0.7 and 7 X 10 3 u.t. Chaotic regions are 
identified with values of the indicator higher than ~ —5 or ~ —4 for the energy surfaces —0.1 and —0.3 or —0.5 and —0.7, respectively. 
See Table [4] for further details. Notice that the values of the indicator are in logarithmic scale. 



torn right panel of Fig. [9] for the energy surface —0.7, the 
portrait shows how the resonances occupy most of the reg- 
ular component. Moreover, major resonances grow within 
the chaotic domains. The most evident resonance is that lo- 
cated at value p XQ ~ 0.5, which has already been identified 
in outer energy surfaces. 

In spite of the short discussion of the p XQ — p ZQ start 
space given above, the efficiency of variational indicators 
to describe the global characteristics of the phase space 
is evident. The VICs easily identify the regular and the 
chaotic domains and their mutual interplay, along with 
phase space mixing processes like resonance overlapping. 
However, the identification of the underlying resonant struc- 
ture by the frequency vectors is not readily understood. We 
should study t he dimensionality o f the tori (e.g. with the 
GALI method. Iskokos et all (|2007h . D12) to identify reso- 
nant orbits, locate the periodic orbits (e.g. with the OFLI, 



iFouchard et al.l ll2002f). D12) and analyse their stability (e.g. 
with the MEGNQ. ICincotta etafl (|2008h ) to look for reso- 
nant orbit families. The process could be time-consuming 
and rather inefficient since the VICs are not the most suit- 
able tools for the task. In fact, the SAMs are far better 
options (see Section [1] and references therein). 

The dynamics of the p XQ —p zo start space can be briefly 
described as a dominant chaotic component based mostly on 
a resonance overlapping regime. This phase mixing process 
produces a more homogeneous chaotic component for the 
energy surface —0.1. 

For further studies concern ing the p xn — p zn start 
space of the ScTS m odel refer to |m 

uzzio et alJ (|2005l ) and 

ICincotta etafl (|200ST ). 

Now, we apply the same package of techniques to the 
xq — zo start space, where the dynamics is totally different 
since it is dominated by the regular orbits. 
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Table 5. For every energy surface in the xo — zo start space, we present the threshold associated to the LI and the estimated percentage 
of chaotic orbits; the threshold associated to the SEiLCE and the estimated percentage of chaotic orbits. 



Energy surface 


Threshold (LI) 


Chaotic orbits (%) -LI 


Threshold (SEiLCE) 


Chaotic orbits (%) -SEZLCE 


-0.1 


1.25 x 10" 4 


~ 3.76% 


9.4 x 10" 6 


~ 4.65% 


-0.3 


5.3 x 10~ 4 


~ 8.15% 


3.2 x 10~ 5 


- 9.2% 


-0.5 


9.4 x 10~ 4 


- 15.06% 


6.5 x 10" 5 


~ 12.74% 


-0.7 


1.7 x 10" 3 


- 16.37% 


9.5 x 10" 5 


~ 16.75% 





Figure 10. Phase space portraits of the SEiLCE for the xo — Zq start space and for four energy surfaces. Top left panel: phase space 
portrait corresponding to the energy surface —0.1 for 998938 i.e. and using a total integration time of 1.17 x 10 u.t. Top right panel: 
-0.3, with 998948 i.e. and 2.4 X 10 4 u.t. Bottom left panel: -0.5, with 999400 i.e. and 1.2 X 10 4 u.t. Bottom right panel: -0.7, with 
998626 i.e. and 7 X 10 3 u.t. Chaotic regions are identified with indicator values higher than ~ —5 or ~ —4 for the energy surfaces —0.1, 
—0.3 and —0.5 or —0.7, respectively. See Table[5]for further details. Notice that the values of the indicator are in logarithmic scale. 



6.2 The xo — Zq start space 

Here, we implement the same procedure applied in Section 
16.11 to study the xo — zo start space. We obtain the phase 
space portraits with the SE/LCE for four energy surfaces, 
calibrating the thresholds with the assistance of the LI, and 
we further discuss the global dynamics of the xo — zo start 
space. 

In order to determine the percentages of the chaotic 
component, we proceed like we did in Section loTTl We present 
the results ordered by energy surface on Table [5] 

In Fig. I lOt we present the phase space portraits of the 
Xo — zo start spaca 11 ! for the same energy surfaces considered 



11 We use a cut-off level around 10~ 10 to have the same scales in 
Fig. 1101 for every energy surface. Then, the percentages of orbits 
that we do not consider are negligible. Those percentages are 



in Section 16.11 Also the total integration times remain the 
same and the VICs used for the study guarantee reliable 
and stable pictures. The numbers of orbits in the samples 
are different and depend on the energy surface: 998938 i.e. 
for the energy surface -0.1, 998948 for -0.3, 999400 for -0.5 
and 998626 for -0.7. 

There is a significant change in the dynamics of the 
Xo—Zo start space with respect to the results for the p XQ —p ZQ 
start space. The corresponding phase space portraits pre- 
sented in Fig. [10] show a dominant regular component (see 
Table [S] columns three and five) instead of the rather large 
chaotic component seen for the p Xo — Pz start space. The 
percentages of the smaller component have a substantial in- 
crement toward the innermost energy surface. Indeed, in this 

~ 0.22% in the energy surface -0.1, - 0.05% in -0.3, ~ 0.02% 
in -0.5 and ~ 0.02% in -0.7. 
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case, the expansion of the chaotic component is due to the 
growth of the different orbital families separatrices. This fact 
is clearly seen in Fig. 1101 On the top left panel of Fig. 1101 we 
can see a very weak separatrix that separates the short-axis 
tubes from the outer long-axis tubes. On the top right panel 
of the same figure, we can see how the former separatrix 
is connected to the others which separate the inner long- 
axis t ubes and the box orbits (see IPapaphilippou fc Laskarl 
( 1998) for further details on the different orbital families and 
their locations in phase space). This becomes evident in the 
bottom panels of Fig. [10] Furthermore, a strong chaotic do- 
main appears in the crossing of the main separatrices. There- 
fore, as the dominant component comprises regular orbits, 
the main contribution to the chaotic component emerges 
from the separatrices of resonances. 

Finally, the dynamics of the xo — zo start space can be 
briefly described as a dominant regular component. The ma- 
jor resonances that separate the main orbital families have 
separatrices which grow when the negative energies increase. 
Furthermore, the growth of the separatrices increases the 
percentage of chaotic orbits very fast. 



7 CONCLUSIONS 

With this report we con c lude a series of investigations (re- 
ported in iMaffione et al.l i|2011al ). Mil and D12) toward an 
efficient choice of a minimal package of techniques to study 
a general Hamiltonian. 

Here, we apply both VIC s and a SAM (the FMFT, 
ISidlichovskv fc Nesvornvl l|l997l )) to study a fairly realistic 
dynamical m o del o f an elliptical galaxy (the ScTS model, 
iMuzzio et all (|2005l )). 

In order to select the appropriate VICs for the ex- 
periments, we not only consider previous comparisons, but 
also ex tend them with new indicators: t he SEZLCE and the 
APLE |Lukes-Gerakopoulos et~aT]|2008l ) . The SEZLCE is an 
efficient estimation of t he 1LCE by means of the MEGNO 
llCincotta fc Sim 3 l2000l ) and works as its reliable alternative 
in case of studying big samples of orbits. The APLE has not 
shown advantages over the FLI/OFLI for the purposes of 
the experiment, but it behaves similarly. 

The SE/LCE seems to improve the performance of the 
MEGNO, in particular situations. Therefore, it should be 
considered in the CIsF presented in D 12 as a suitable alter- 
native. 

On the other hand, we consider the performance of the 
FMFT as the representative SAM to compare with the VICs. 
The F MFT is an im provement of the FMA outlined by 
Laskar (|Laskarlll990l ). which is widely used by the scientific 
community. 

The main advantage of the SAMs is a fast computation 
of the frequencies. However, the performance of the FMFT 
as a global indicator of chaos is not as efficient as the perfor- 
mances shown by the VICs in the experiments. The SAMs 
are designed to compute the frequencies, which are quanti- 
ties strictly related to regular motion. Using such techniques 
as global indicators of chaos involves forcing the methods to 
do something for which they are not designed. The VICs are 
based on the concept of local exponential divergence. Hence, 
the detection of chaos is their main purpose. A natural con- 
sequence of this is that an efficient application of those tech- 



niques is based on their complementary implementation (see 
Section [1] for references) . The SAMs are of use to describe 
the resonance web while the VICs to study the interplay of 
regular and chaotic domains (see Section [6}. 

Finally, the recommended CIsF for the analysis of a gen- 
eral Hamiltonian is composed of the MEGNO/SE/LCE, the 
FLI/OFLI and the GALI 2 jv: the MEGNO/SEZLCE and the 
FLI/OFLI as global indicators of chaos to obtain the phase 
space portraits and display the interaction between regular 
and chaotic orbits; the GALI2JV to analyse small samples 
within regions of complex dynamic s and regions wher e the 
chaotic component is dominant fsee lSkokos et all |2007l ) and 
D12). Furthermore, the CIsF can be used with a SAM (e.g. 
the FMFT) in order to characterise the resonance web with 
the frequency vectors. 
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